library(DESeq2) 
library(ggplot2)  
library(ggrepel) 
library(here)
library(dplyr)
library(plotly)

input <- "mRNA_deseq_results.RData"
load(here("processed_data", input))

PCA

Corrected PCA

Raw PCA

3D PCA

vsd <- vst(dds) # variance stabilizing transformation of Deseq object
pcaResult <- prcomp(t(assay(vsd))) # transpose vsd and perform PCA
pcaData3D <- as.data.frame(pcaResult$x[, 1:3]) # Extract the first three principal components
pcaData3D$condition <- vsd$condition # New condition column is added

# Calculate the percent variance for the first three components
percentVar <- round(100 * pcaResult$sdev^2 / sum(pcaResult$sdev^2)) 

# Create a 3D plot with Plotly
fig <- plot_ly(data = pcaData3D, x = ~PC1, y = ~PC2, z = ~PC3, color = ~condition, 
               type = 'scatter3d', mode = 'markers', text = rownames(pcaData3D))

fig <- fig %>% layout(
  title = "3D PCA Plot",
  scene = list(
    xaxis = list(title = paste0("PC1: ", percentVar[1], "% variance")),
    yaxis = list(title = paste0("PC2: ", percentVar[2], "% variance")),
    zaxis = list(title = paste0("PC3: ", percentVar[3], "% variance"))
  )
)

fig